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Abstract 

Quasi-one-dimensional spin-Peierls and spin-ladder systems are character- 
ized by a gap in the spin-excitation spectrum, which can be modeled at low 
energies by that of Dirac fermions with a mass. In the presence of disorder 
these systems can still be described by a Dirac fermion model, but with a ran- 
dom mass. Some peculiar properties, like the Dyson singularity in the density 
of states, are well known and attributed to creation of low-energy states due to 
the disorder. We take one step further and study single-particle correlations 
by means of Berezinskii's diagram technique. We find that, at low energy e, 
the single— particle Green function decays in real space likeG(x,e) oc (1/x) 3 / 2 . 
It follows that at these energies the correlations in the disordered system are 
strong - even stronger than in the pure system without the gap. Addition- 
ally, we study the effects of boundaries on the local density of states. We 
find that the latter is logarithmically (in the energy) enhanced close to the 
boundary. This enhancement decays into the bulk as and the density 

of states saturates to its bulk value on the scale L e oc ln 2 (l/e). This scale is 
different from the Thouless localization length A e oc ln(l/e). We also discuss 
some implications of these results for the spin systems and their relation to 
the investigations based on real-space renormalization group approach. 

I. INTRODUCTION 

One-dimensional quantum electron (and spin) systems have attracted considerable at- 
tention of theorists over the decades. The interest to these systems have been particularly 
boosted in recent years by the rapid development of experimental techniques. The latter in- 
clude the discovery of various non-organic quasi-one-dimensional compounds. In particular, 
the materials we shall be concerned with in this paper are the recently discovered GeCuO^ 
spin-Peierls compounds^ and the spin-ladder compounds (^0)2-P207 and SrCu 2 0$. The 
modern experimental techniques allow measurements on these inorganic compounds, which 
were either impossible or inaccurate with the organic spin-Peierls materials in the past. 
These measurements not only can be now performed on single crystals, but also involve 
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controlled doping by impurities. The latter possibility brings up the old issue of the disorder 
effects in one-dimensional systems: a fascinating subject which was considered in the past 
as being somewhat academic. 

Among various one-dimensional disordered systems, random exchange spin chains have 
been studied, using the so-called real-space renormalization group, with very intriguing 
resultsIF! In particular, it has been shown! that the typical behavior of the correlation 
functions may be quite different from the average behavior, which is more relevant from 
the experimental point of view. For instance, while the typical correlation functions usually 
decay exponentially, the average ones can be power law decaying. Moreover, even when 
the average correlation functions decay exponentially, the correlation length is different 
(bigger) from that of the typical correlationsB. This feature is a consequence of rare disorder 
configurations dominating the long distance behavior of the correlation functions. The 
analogy of these spin systems with the effective disordered fermionic models have not been 
investigated in detail. 

In this paper we shall partly fill this gap, and investigate more deeply the fermionic 
models. The main reason is that fermions naturally appear in the description of pure spin- 
Peierls chains and especially in the spin ladders. Moreover, this alternative approach is 
somewhat complementary to the real space renormalization group analysis, thus together 
they would provide a complete and satisfying description of such disordered systems. 

Specifically, considering that both the spin-Peierls and the spin-ladder systems have a 
spin gap in the excitation spectrum, it is not surprising that the effective fermionic model for 
both systems is a model of massive Dirac fermions. The effect of non-magnetic impurities is 
to randomize the mass (see also the next Section). Thus, the ultimate fermionic model, on 
which we focus in what follows, is the one-dimensional random mass Dirac fermion model. 

A great deal is known about simple self-averaging quantities for this model, like the total 
density of states and the localization length. 

For a Gaussian (white noise) distribution of the mass variable m(x) the density of states 
was calculated by Ovchinnikov and Erikhmani. Using Fokker-Planck type equations for the 
probability distribution of the wave-function phase, they obtained the divergent density of 



as the energy approaches the midgap e — > 0. 

Physically, the appearance of the singularity in the density of states can be easily un- 
derstood. Indeed, for a single kink (sign change) of the mass, there always is a zero-energy 
bound state, localized around the point in space where m(x) = 0. For many kinks, there 
are many such localized states (also referred to in the literature as zero-modes or solitons). 
If they were not overlapping, the density of states would have a 5-function peak at e = 0. 
In fact the zero-modes are overlapping. Hence the 5-peak broadens, but the singularity re- 
mains - Eq. ([!]). [In order to give a qualitative explanation as to why the broadened 5-peak 
takes exactly the shape (P, a more elaborated argumentation is needed: see, e.g., Ref. |7|.] It 
must be noted that the singularity in the density of states of the form (|l|) has been discov- 
ered by Dyson back in 19530 for a model of a disordered harmonic chain. In the electronic 
spectrum at the center of the Brillouin zone such a singularity was identified by Weissman 
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and Cohan! for the case of an off-diagonal disorder (random hopping model). The latter 
model is in fact directly related to the random mass Dirac problem. 

According to the Thouless relation^, the localization length corresponding to the density 
of states (H|) takes the form 

A e ~nn(l/e), (2) 

where I is the mean free path. The fact that the localization length diverges in the middle 
of the gap makes the random mass Dirac problem different from the usual one-dimensional 
localization problem. For the latter, the localization length essentially coincides with the 
mean free path I, so that there is only one length scale in the problem. Consequently, it 
makes sense to consider only the distances x ^> / where the correlation functions decay 
exponentially^. On the contrary, for the random mass Dirac fermions, various length scales 
come into play: the mean free path I, the localization length A e , and what we shall call the 
correlation length L € . We are going to argue that 

L e ~ nn 2 (l/e) . 

(This is in agreement with the renormalization group studies, mentioned above; see also 



Section |VTl) The correlation length diverges even faster than the localization length. While 
we do expect that the correlation functions for the random mass Dirac fermion model behave 
in the same way as those for standard one-dimensional disordered systems at distances 
x ^> L 6 , there is a different regime I < i < L £ , arising at low energies. A new physics 
emerges in this regime, the understanding of which this paper is intended to contribute to. 
Specifically, the system seems to exhibit some kind of criticality - all the known correlation 
functions are power-law decaying with universal exponents. This behavior is determined by 
rare fluctuations when the particle wave-functions have large amplitudes. 

As we already mentioned, the simple self-averaging quantities, like the total density of 
states, can be computed by various methods like the Fokker-Planck type equations for the 
probability distributions or the replica tricklH. All these methods are, however, ill-suited 
for determining the correlation functions. In their recent preprint Shelton and Tsvelik 
succeeded to make use of an elegant mapping of the Dirac problem at zero energy (e = 0) 
onto the so-called Liouville quantum mechanics^ for calculating the disorder averages of the 
products iPq(xi)iPq(x2)...iPq(xn) of the zero-energy wave-function amplitudes. While this 
mapping does provide an interesting insight into the random mass Dirac problem (strictly) 
at zero energy, we are interested to get more information about the correlation functions, in 
particular about the correlations involving the phase of the wave-functions and about their 
energy dependence. 

On the other hand, a method, which does (in principle) allow to calculate all the correla- 
tion functions, has been around for almost 25 years. The method is the Berezinskii diagram 
technique originally invented (back in 1973) to confirm the Mott hypothesis of a; 2 In a; vanish- 
ing conductivity for the standard one-dimensional localization problem^. The Berezinskii 
method was extended by A. A. Gogolin and Mel'nikov (GM)0 to the case of a half-filled 
electron band with a random backscattering, the model equivalent to the random mass 
Dirac fermions. Needless to say that GM were able to reproduce the Dyson singularity 
in the density of states. More importantly, they calculated the current-current and the 
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density-density correlation functions. The current-current correlation function has no ap- 
parent physical meaning for the spin systems. From the point of view of the latter, the 
density-density correlation function is more interesting, for it turns out to be proportional 
to the staggered susceptibility (at least for spin-Peierls systems, see also below) . GM found 
that, in the region ICs, the contribution to the density-density correlation function at low 
frequency u of the state with energy e such that i<L e decays as 

1\ 3/2 



(^U(x)^(x)^(0)^(0)) cc (-J . (3) 



(This asymptotics was reproduced by Shelton and TsvelikS for e = uj = 0). 

In this paper we further study the random mass Dirac problem following the line pi- 
oneered by GM. After defining the model and discussing its relation to spin systems, we 
briefly introduce the Berezinskii technique (we feel that that would be convenient for the 
reader since a transparent explanation of the Berezinskii technique is often neglected in the 
literature). Then we apply the technique to the calculation of the single-particle Green 
function at low energy. Amusingly, we find that the Green function decays as (1/x) 3 ^ 2 (for 
I < a; « L e ). Namely, it is characterized by the same exponent as the density-density 
correlation function. 

In the second part of the paper we study a somewhat different (though related) question. 
Namely, we calculate the local density of states in the presence of a boundary. Unlike the case 
of the usual localization [for which the effects of boundaries were studied by Al'tshuler and 
Prigodin (APft for the random mass Dirac problem one might expect that the boundaries 
affect already the average local density of states (and not just its distribution function). 
Indeed, it turns out that the Dyson singularity is enhanced close to the boundary as ln(l/e), 
as compared to its bulk value. The local density of states then decays into the bulk of the 
sample following the law l/\/x, and finally saturates to the bulk value at x ~ L t . We also 
confirm these results by a real-space renormalization group analysis. Some consequences 
for the thermodynamics of the spin systems are discussed in the last Section. 

We would like to stress that, despite all the results mentioned above, we are just starting 
to study the problem in this paper. Ultimately, a comprehensive understanding of the 
new physics, arising in the random mass Dirac problem will only be achieved when the 
fluctuations of the correlation functions are understood. The problem is mathematically 
involved, so that, unlike the usual case@, even the distribution function of the simplest 
quantity, the local density of states, is unknown at present. We hope to return to these 
questions in future publications^. 

II. THE MODEL 

The one-dimensional Dirac Hamiltonian is of the form 

H D = J dx{v^(-id x )a z y + m(x)^a x ^}, (4) 

where a x , a z are the Pauli matrices, and \E^(x) = R^(x), L\x) . The operators R(x) [L(x)] 
stand for the chiral right (left) moving electron fields. The mass term m(x) is randomly 
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distributed. Throughout this paper we shall assume m(x) be Gaussian with vanishing 
average: (m(x)) = 0. Its distribution is therefore determined by the correlation function 

v 2 

(m(x)m(x')} = ~$( x ~ x ') i (5) 

where I is the mean free path (and v is the velocity). 

In practice, the Dirac model (Q) usually arises as a continuum version of a lattice model, 
so that R and L are the chiral components of the continuum limit for a lattice electron field 
operator 

i) n -> e ikFX R(x) + e- ikpx L(x), for x = na . (6) 

Here ao is the lattice spacing and the condition of being at the half-filling for the lattice 
model means that 4/cpao = 2tt. From (|6]) it is clear that (|j) is just another way to write 
down the random backscattering electron model at half-filling. 

Notice that there is no forward scattering in the problem. This is rigorously true for a 
simple tight-binding model at half-filling with random nearest-neighbor hopping integrals. 
This model is particle-hole symmetric 

which prevents any generation of forward scattering processes. Hence this model in the 
continuum limit would reproduce the Hamiltonian (|4]). 

The Hamiltonian (f|) has been shown to also describe the low-energy limit of one- 
dimensional spin systems with a spin gap in the excitation spectrum (see Refs [T7|,[n|). For 
the convenience of the reader we shall present the qualitative side of the argument. 

• Spin-Peierls systems. The spin part of the Hamiltonian for a single spin-1/2 chain 
interacting with phonons takes the form 

HsP = X! JSnSn+l + X! ^n( — ~L) n Sn.S n +l , (7) 
n n 

where S n are spin-1/2 operators, J is the exchange coupling constant, and A n measures 
the strength of the dimerization caused by the interaction with the phonons. Below 
the spin-Peierls transition temperature, one usually assumes A n — ► (A n ) = A (the 
mean-field approach is appropriate for the phonons normally have a three-dimensional 
dispersion relation). The passage to the model (f|) is as follows. The XY-version of 
the model (0) can be mapped onto a fermionic model in a standard way by using 
the Jordan-Wigner transformation. Taking the continuum limit leads to (§) with the 
parameter A determining the mass m . The J z coupling corresponds to a four-fermion 
interaction term, but it is irrelevant because of the spin gap. The main effect of the 
doping is to randomize the mass term m(x). As argued in Ref. |T7|, a single impurity 
effectively introduces a domain wall into the system - a kink separating phases with 
different signs of the dimerization order parameter. The mass variable, proportional 
to the dimerization strength, will thus also have a kink at the domain wall position 
Xq, changing from m{x) = mo for x < Xq to m(x) = —mo for x > Xq. 
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Spin-ladder systems. The Hamiltonian of two coupled spin-1/2 chains is 

Hsh = JSn,jSn+l,j + ^ J±S„.,lSn,2 , (8) 

n,j=l,2 n 



where the coupling between the chains, J±, is assumed to be small (this it not quali- 
tatively restrictive but allows one to construct a consistent theory). In the continuum 
limit description, the spin density operators on each chain have slow-varying compo- 
nents as well as fast (staggered) components: 

S nd -> S(x) + (-l) n n,(x) . (9) 



The staggered correlations are stronger than the uniform correlations, so the most 
relevant interaction between the chains is 

J '±clq 1 J dxn\(x) • n^x) (10) 



(other couplings are, in fact, marginal). The scaling dimension of the coupling (TO) is 
equal to unity. Therefore, as it was realized by Shelton, Nersesyan, and TsveliktB, this 
coupling can be conveniently re-fermionized. The result is the Dirac model (f|) with 
the mass parameter proportional to coupling J^0. Let us now discuss what the effect 
of doping would be. Consider, for example, the situation in LaiCu0 2 .5 spin-ladder 
system doped by Sr (which substitutes La)c3. For low Sr concentrations (x < 0.02), 
the holes (carrying spin-l/2's) are localized on the Oxygen atoms, i.e. in between 
neighboring sites of the Cu magnetic lattice. As a result the magnetic sites from the 
right of such an impurity are effectively re-numbered (shifted by one). Consequently, 
the staggered magnetization changes sign causing a sign change of the coupling fliCf ) 
and of the mass term of the corresponding Dirac problem!!. 

Thus, the above considerations of the doping effects lead, for both the spin-Peierls and 
the spin-ladder systems, to the random mass Dirac model. Specifically, the mass turns out 
to be distributed as a random 'telegraph signal' 

m(x) = m sgn(x — X{) , (11) 

i 

X{ being randomly distributed. At low energies (namely at energies well below the gap: 
e rn ), the model is still equivalent to the one with a Gaussian distribution of the massQ. 
Since we are mainly interested in the asymptotic behavior of the correlation functions (i.e., 
in the parameter region where the correlation functions are expected to be universal: e tuq 
while I C i < L e ), we can safely assume the Gaussian distribution. Then the correlation 
functions can be calculated by means of the Berezinskii diagram technique0 introduced in 
the following Section. 
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III. BEREZINSKII DIAGRAM TECHNIQUE 



The Berezinskii technique is based on the perturbation theory for the one-particle 
Schrodinger equation (rather than on the many-body perturbation theory). Consider a 
translationally invariant one-particle Hamiltonian H perturbed by a random potential (ran- 
dom mass) term rh. The bare retarded Green function in the energy representation is given 
by 

G 0+ = i . (12) 

e-H + i5 

The exact Green function for a given realization of the disorder is defined as 

G + = — ^ . (13) 

e — Ho — rh + i8 

This Green function satisfies the equation 

G + = G 0+ +G 0+ mG + . (14) 



Iterating (|T4| ) one obtains a standard perturbative expansion for the Green function. A 



typical term of the perturbative series is of the form 

G 0+ mG 0+ mG 0+ ...G 0+ mG 0+ . (15) 

The next step is to take the disorder average. Clearly only the 2n-th order terms survive 
(i.e., the terms containing an even number of impurity potentials). The disorder average 
(m(xi)m(x2)--.m(x2n)) is then given by a sum of products of all possible pair correlators (|J). 
(An additional matrix structure of the mass term simply ensures that the left- and right 
movers are interchanged in every scattering event.) For example, at second order there is 
just one term 

/ cfeiG 0+ (x\xi)G 0+ (xi,xi)G 0+ (xi,x) (16) 

(The usage of the coordinate rather than the momentum representation is vital, as it will 
become clear shortly.) At fourth order there are two terms and two integration variables, 
etc.. 

Let us now represent these contributions to the Green function graphically. The bare 
electron Green functions are represented by solid lines. The disorder potential correlators 
by 'interaction' (wavy) lines, which shall always be vertically aligned. The term flTB] ) for 
x' < Xi < x is then given by the graph Fig.|l|(a). This graph is also an example of an 
elementary 'interaction' vertex. Other contributions are represented by more complicated 
diagrams where several vertices Fig.|I](a) and (or) other vertices are joined. An elementary 
consideration of the scattering processes described by the model (|]) shows that Fig.|T|(a-d) 
exhausts all the possible vertices. Using the many-body terminology, there are no forward 
scattering vertices [simply because there is no scattering within the same chiral branch in 
(||)]. The vertices Fig.|I](c,d) are referred to in the literature as umklapp scattering vertices. 
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It is very important that, by construction, the diagrams are x-ordered: X\ < ...X{ < x' < 

Xi^-\ <C ... <C Xj <C X <C <C ... <C -^2n^^* 

As we have already stressed, each scattering event interchanges the right- and the left- 
moving electrons. The bare Green functions for them are given by 



G°^{x', x- e) = - % -Q(x' - x)e l ^ x ~ xl) , G° L + {x', x; e) = - % -d(x' - x)e l ^ x '~ x) , (17) 



where 8(x) is the step function. Alternatively, one can use the bare Green function of the 
combined electron field 

G 0+ (x', x; e) = e ikF{x - x,) G ^{x', x; e) + e^^^G"^', x; e) = -^ e l(kF+ ^ x - x ' 1 . (18) 

The next important observation is that, since all the diagrams are space-ordered, the 
Green functions (OJTSp are actually factorized in each given graph. This enables one to split 
the unperturbed Green function G 0+ (xi, xf, e) = — ^exp i{kp + — Xi) (for Xj > Xi) 

into two factors exp — i(A;^ + ^)x{ and exp +i(kp + ^)xj , and to assign the coordinate 
dependence only to the vertices and to the end points, not to the lines. As a result, the solid 
lines just tell us how the vertices are joined: the analytic expression for any graph is given 
by a product of the vertex (and end point) factors. Note that the factors corresponding to 
the vertices Figfj] (a,b), -1/(2Z)@ and —l/l respectively, are ^-independent. The factors 
(— l//)exp [±i(e/v)x] corresponding to the vertices Fig.[TJ (c,d) also weakly depend on x (as 
e — > 0). The reason is that the Fermi wave-vector kp corresponds to half-filling: e tAkFX = 1. 
As it was realized by GM0, these (umklapp) vertices make the problem different from the 
standard localization problem. (At large energies, i.e. well away from half-filling, these 
vertices strongly oscillate and give a negligible contribution after integration over the spatial 
coordinates is performedS.) 

So far we have merely reformulated the problem. However, the cardinal idea of the 
Berezinskii methods is to provide a way of summing all the diagrams. The procedure is 
based on an ingenious classification of diagram Since each diagram is characterized by a 
fixed number of vertices in a definite spatial order, one can split it into a product of three 
integrals over the vertex contributions: the contributions from the region to the left of x', 
those made in the region to the right of x, and those made in the central region (between 
x' and x). These integrals, which are referred to in the literature as 'Berezinskii blocks', 
are then classified by the number of lines passing through their boundaries. For instance, 
R m (x') is the sum of all contributions to the left of x' with 4m lines, and, equivalently, 
R m {x) defines the sum of all contributions to the right of x. By joining interaction vertices 
to Berezinskii's block of a given order, one finds Berezinskii's blocks of higher (and lower) 
order. This allows for recurrence relations to be derived. Instead of trying to give here a 
general description of the procedure, we shall illustrate the method at work while calculating 
the single-particle Green function in the next Section. 



IV. SINGLE PARTICLE GREEN FUNCTION 

By multiplying the iteration series (|i~5|) with different energies (and spatial coordinates) 
and applying the Berezinskii method one can, in principle, calculate all the correlation 
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functions. So, the density-density (and current-current} correlation functions, which are 
two-particle correlation functions, were found by GMEj (see Section VI for discussion). 
In this paper we concentrate on the single-particle properties of the random mass Dirac 
problem. We thus apply the Berezinskii technique for calculating the single particle Green 
function G + (x f , x; e). 



A. Basic equations 

We start classifying the graphs contributing to the single-particle Green function by 
considering all possible end point configurations. Since the electron line can be drawn 
either to the left or to the right of each of the two end points x', x, there are four possible 
configurations, as shown in Fig.^. In going along a diagram, the number of the electron lines 
passing through the current vertical cross section of the digram may be changed (by some 
amount Ag) when one encounters an interaction vertex. So, the vertices Fig.|l|(a,b) do not 
change the number of lines, Ag = 0. The vertices Fig.|l](c,d) do change the number of lines by 
Ag = ±4. An elementary analysis shows that the end point configurations Fig.^(c,d) would 
require Ag = ±2. This is impossible with the set of elementary vertices at our disposal. So, 
the configurations Fig.^|(c,d) do not contribute to the Green function. 

We are therefore left with the two end point configurations Fig.^|(a,b). The structure of 
the Berezinskii blocks corresponding to these configurations is shown in Fig.|3|(a,b). The R 
and R blocks can only have multiple of 4 number of the electron lines. They are therefore 
identical for both configurations Fig.^|(a,b). On the other hand, the central blocks Z^ a > and 
differ. Indeed, given the left-hand-side block R m i has 4m' lines and the right-hand-side 
block R m has 4m lines, the central block Z a must have 4m' + 1 (4m + 1) lines on the left 
(right) while the central block Z h must have 4m' — 1 (4m — 1) lines on the left (right). [Thus, 
m > for a-configuration while m > 1 for 6-configuration.] 

In terms of the Berezinskii blocks the single-particle Green function is expressed as 
follows (x' < x) 

G+(x',x)= --{ £ R m ,(x')Z { J m (x',x)R m (x) 
v I — 

^ m,m'=0 

oo I 

+ ]T R m/ (x')Z^ !rn (x' 1 x)R m (x)\ . (19) 

m,m'=l J 

Consider now the Green function for coinciding points x' = x — 0. When the end points 
merge, the region of the spatial integration in the central blocks collapses so that only the 
trivial diagrams for Z^'^ survive (those which do not include any interaction vertices). 
Hence the boundary conditions 

Z { ^ b l(x -0,x) = <W, m for m > 0(1) . (20) 

For an infinite sample, the boundary conditions to be imposed on i? m 's are 

R (x) = R (x) = 1 . (21) 

Thus, the equal-point Green function is given by 
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G + (x-0,x) = — \l+2 y £R m (x)R m (x) \ . (22) 

V { m=l J 

The local density of states is denned as 

p(x,e) = --lmG + (x - 0,x;e) . (23) 

7T 



In an infinite system the quantities ([22]) and (|23|) do not actually depend on sc, for the 
disorder average effectively restores translation invariance. (The situation is different for 
systems with boundaries we discuss in Section |V|.) Note the appearance of the factor of two 
in front of the sum in the equation (^2]). This factor is due to configurations Fig.Q mirror 
imaging each other when the end points merge, as illustrated in Fig.[5[ 

The expression (|23|) was used by GM for calculating the density of states. They found 
the following differential recurrence relation for the Berezinskii blocks R m 



dR m {x) 

dx I 



- {4m 2 R m {x) + m(2m - l)e~^*£ r , 
+ m(2m + l^-^Rm+^x)} . (24) 



The right-hand-side blocks R m {x) mirror the blocks R m (x), so they satisfy ( |24|) with 
x — > —x. An integral equation equivalent to the equation (|24]) is obtained by joining the 
elementary vertices to R m (x) in all possible ways, but avoiding creation of closed electron 
loops. (The loop diagrams are not allowed because we study a one-particle problem. Put 
another way, they vanish because they involve a retarded Green function.) The process is 
illustrated in Fig.f|. The first term in (|24]) comes from joining the vertices Fig.|I](a,b) 

4m 2 = ^Am + 2m(2m - 1) . 

The second and the third term derive from joining the vertices Fig.[I](c) and (d) respectively. 

Carrying out an analogous calculation (Fig.fD for the central blocks Z^ a ' b \ we obtained 
the differential recurrence relations satisfied by these blocks (in computing the integer coef- 
ficients one must be careful as to not allowing for the electron loops in the whole diagram): 

^-^•iW-^^^^W (25) 

and 

= - 8m '- a 4m+1 Zg,(*) (26) 

- ("-i)(2"-i) e -,^g,_ iM _ a^si+Ji^zss. «m 

Only the active spatial variable x is shown (the variable x' is suppressed), * stands for the 
index w! which plays a role of a parameter. Defining the Berezinskii blocks [c.f. (|i~9|)1, we 
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have incorporated the end point exponentials exp (^Li^x^j into Z^ a,b \ Hence the first terms 
on the right-hand-side of the above equations. [The absence of e ±lkFX factors means that we 
are actually calculating Gr, the total Green function can be restored by using the formula 

©•] 

As a next step it is instructive to pass to the momentum representation for the Green 
function. Before doing so we notice that the spatial dependence of the side blocks can easily 
be separate dB. Indeed, the substitution 

R m (x) = (~l) m e- l ^ mx R m (27) 

reduces the equations (pip and (|2T| ) to the algebraic recurrence relation 

isR m = AmR m - (2m - l)i? m _i - (2m + l)R m +i (28) 

for the x-independent quantities R m obeying the boundary condition 

R = 1 . (29) 

We have introduced the dimensionless energy variable 

Ael . . 

s = — . 30 

v 

It will also be of convenience to pass to the dimensionless spatial coordinate 

» = f. (3D 

Now we define the objects 
Q%»{k,s) = l(-l) m £ (-1)™' / dye-^y-y'^ m y~ m ^zl:; b l(ly',ly)R mf . (32) 

m'=0(l) J y , 

The sum over m' in the above formula starts with for a-configurations and with 1 for 
^-configurations respectively. 

In terms of the objects (|3"2T), the Fourier transformed Green function takes the simple 
form 



G R (K,s)=lJ dye-^y-y^G R (ly'Jy;s) 



v 

oo 



£ R m Q$(K,s)+ £ R m Q$(K,s) 



_m=0 m=l 



(33) 



Finally, applying the operation ([32]) to Eqs. (|25|j26| ) and using the boundary conditions 
), one can derive algebraic recurrence relations for Qm' i K i s )- We found 

is (m + \) Qt ] - WW - 8m2+ 2 4m + 1 QL Q) + 

m(2m - l)Q ( Si + I) (2m + l)Qtli = R m (34) 
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and 



• ( l \n(b) ■ nib) 8m 2 -4m + l , b) 
tsl^m--) - ikQ { J Q { m > + 

(m - l)(2m - l)Q<*-i + m(2m + = i? m , (35) 

where we have suppressed the inactive variables (k, s). 

B. Asymptotic solution of the equations 

In order to compute the Green function we need to solve the recurrence relations 
(|8|,[34],[35|) for the quantities R m and Q^' b ' and use the formula (|33|). 



As we already mentioned, the recurrence relation for -R m 's was derived by GM. 
They also discovered that the generating function R(C) = J2 m Rm( m satisfies a second order 
differential equation that happens to be of a hypergeometric type. It can therefore be 
explicitly solved by, e.g., contour integrals method. The solution found by GM rea ds0 

Rm = K Q {-is/A) Jo jpffci (36) 

where K is the MacDonald function. (It can be checked that the boundary condition Rq = 1 
is fulfilled.) 

The recurrence relations (I31J5B) for Q^'^'s are of a more complicated nature. So, the sec- 
ond order differential equations satisfied by the corresponding generating functions Q^ a ' b '(C) 



are more complicated than the hypergeometric equation. (These equations have four singular 
points; we have derived them but we do not feel that they are in situ here.) Unfortunately, 
we have not been able to solve these equations. 

On the other hand, we are mostly interested in the low-energy regime s — > 0. There one 
can simplify the recurrence relations and change over to differential equations. Indeed, in 
the low-energy limit the high order diagrams matter, so that we shall restrict ourselves to 
the parameter region 

s — ► 0, m — > oo while sm finite. (37) 
An inspection of the formula (|36|), shows that in the limit fl37j) 

*« = -SF5) Mz) • (38) 

where we have introduced the variable 

z = V ' —2ims . (39) 

[For large m we can replace [t/(t + l)] m by exp(— m/t) in (^6|) and the remaining integral, 
upon rescaling the integration variable, gives the MacDonald function.] 

An alternative way to obtain this approximate solution is to pass to the limit of large m 
in (|28|) B. Expanding R m ±i and neglecting the terms vanishing in the limit (|37|), one finds 
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m , 5 + m— h -smR m = 0. 

dm z dm 2 

Making the substitution fl3D|) one arrives at the modified Bessel equation 

(PR ldR _ R _ 

d 2 z z dz 

A general solution to this equation is R(z) = AKq(z) + BIq(z). The coefficient B must 
vanish for the function Iq(z) diverges at infinity. The coefficient A is to be determined by 
the boundary condition fl29|) . Since the MacDonald function diverges as z — > 0, we must 
recall that the variable z is discrete. There is a z min = \/—2is, so that A = 1/Kq{z v ^) ~ 
— 2/ln(— is), in agreement with d38[) within the leading logarithmic accuracy. (Within this 
approximation one must generally neglect C-numbers as compared to the log's. Yet we are 
keeping the imaginary part of the log's. Alternatively, the imaginary parts can be restored 
at the end of the calculation by a simple analytic continuation, i.e. demanding that the 
Green function be retarded.) 

Let us now determine the quantities Q^'^ in the limit (j37D . The recurrence relations 
(|34|) and fl3"5| ) become identical for large m and reduce to 

2m ~^hr + 4m^^ + ( i - iK + ism) = R m (40) 

dm 1 dm \2 J 

Since the summation over m in formula (|33D is now to be changed to the integration over 
z [according to J2m=o(i) ~ * (V s ) fo° z dz], the Green function is determined by the sum 
Q = Q( a ) + Q( fe ). The equation for the latter follows from (||) and (0): 

z2 d*Q + ^dQ _ 2 + ^ _ = 
dz z dz 



Rewriting the equation ([y]) in the form 

d{zQ) 



d 

dz 



z- 



dz 



9i 'K \ 

z+—j(zQ)=AR(z), (42) 



we identify it as an inhomogeneous modified Bessel equation for the function zQ(z). 

Equations of this type are usually solved by means of the Lebedev-Kontorovich 
transformation^. We employ the Lebedev-Kontorovich transformation in the form 

oo 

Q(r) = J dzK lT (z)Q(z) , (43) 
o 

where K ir (z) is the MacDonald function of a purely imaginary index. The transformation 
inverse to flUf) is 

oo 

Q{z) = 4- / drr sinh^r) K tT (z)Q{r) . (44) 
n z z J 
o 
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Applying the integral operator ([43]) to the equation (f42|), re-grouping the terms and then 
using the inverse transformation ([44"1), we find the solution in the form 



Q(«,z) = ^/ rfrrsinh(vrr)-^^- J d^ K lT (O^R(0 , (45) 



r 2 + 2in 




where we have restored the dependence of the function Q on the parameter k. 
At this point we notice that (|45|) is normalized as 



d,K 



J ^Q(k,z) = 2R(z). (46) 

—00 

According to the definition fl3~2"|) , the normalization (^) corresponds to the boundary con- 
ditions (gD|). 

Substituting the asymptotic form fl38|) into ( ]45| ) and looking the ^-integral up in Ref. |27], 
we finally obtain 

OO 

Q(k, = frfrr tanh f —) J^L . (47) 

v ' 



C. Calculation of Green's function 



It is now instructive to Fourier transform the Green function ( p3| ) back to the coordinate 
space: 

G + R {y,s) = - dzzR(z)Q(y,z) . (48) 
vs Jo 

The function Q in the space domain is of the form 



00 



Q(y,z) = - -^-— jdrrtanh K iT (z) e -^ . (49) 







The z-integration is performed in Appendix A. The final result for the Green function 

is 

G%(y,a)= ; t ■ M y), (50) 

vsm {—is) 

where the function F(y) is defined by 

Z" 00 2 2 

F(y) = 2 y da; a; tanh a; (1 — tanh 2 x)e~^ lVX . (51) 
The limiting behavior of the function F(y) is 
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« > f 1 - \y + off 2 ) ■ for » « 1 

F <»> = {^( i )« a + oto*»),fa B »i < 52 > 

The function -F(y) is plotted in Fig.|| Notice that the correct density of states (Dyson 
singularity) is reproduced as y — > 0: 

P(s) = --ImGUv - 0, a) * .%Z S ' ( 53 ) 
7r sin {l/sj 

where po — i/(2irv) is the bare density of states (in the chiral branches of the fermion 
spectrum) . 



V. BOUNDARY EFFECTS 

In this section we study the influence of boundaries onto single-particle properties of 
the random mass Dirac problem. This study is largely inspired by the investigation of the 
boundary effects for the standard one-dimensional localization by AP succeeded in 

calculating the whole distribution function of the local density of states (for an infinite 
system as well as in the presence of boundaries). Technically, this success was due to a 
factorization of the high-order correlation functions into combinations of known Beresinskii 
blockHS. Unfortunately, for the random mass Dirac problem such a factorization does not 
occur. So, at the time being we do not know the distribution functions. On the other hand, 
as we shall see shortly, already the simplest quantity - the local density of states - displays 
quite an interesting boundary behavior for the Dirac problem. 



A. Modification of the diagram technique 

We consider a semi-infinite sample (x > 0), i.e. we assume an infinite potential wall to 
be situated at the origin (the electrons can not leave the sample). As a result of the total 
reflection from the boundary the single-particle eigenf unctions of the Hamiltonian without 
the disorder potential change over from the plane waves to the standing waves: 

1 • [2 

Q i P x _ > i / _ sin(px) . 



y/L V L 

The bare Green function therefore takes the form 

1 r 



G°+ en (x', x;e) = - G 0+ (x - x', e) + G 0+ (x' - x, e) 



- G 0+ (x + x', e) - G 0+ (-x - x', ej 

= G 0+ (x-x',e)-G 0+ (x + x',e) . (54) 

The last equality follows from the fact that G 0+ (— x) = G 0+ (x) for the translationally 
invariant system, Eq.([] 



Let us take a closer look at the second term in the expression (54) which is due to the 
presence of the boundary. This term can be written as 
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- G 0+ (x + x',e) = y^+*M(*+*') = {_!} |yCle i ^+ e /°) a: '| |yCl e l '(^+ e / v )4 . (55) 

In terms of the Berezinskii diagrams this term can be represented as shown in Fig.[7]. The 
second and the third factors in curly brackets in the expression (|55"D correspond to the 
external vertices x' and x respectively, while the {—1} factor must be attached to the new 
element of the diagram technique - the 'boundary vertex' at x = 0. Physically this vertex can 
be interpreted as a total reflection from the boundary. (Any phase factor, like — 1 = exp(i7r), 
corresponds to a shift of the position of the boundary and can be omitted.) 

Thus, in the presence of a totally reflective boundary the Berezinskii diagram technique 
is modified by adding the boundary vertex: Fig.|7|. 



B. Density of states 

The density of states at the point x is given by the general formula (0). We must 
therefore consider the one-loop diagrams: the electron line starting at the point x returns, 
after undergoing scattering processes and reflections from the boundary, to the same point x. 
An example of a diagram involving the boundary scattering is shown in Fig|| Notice that 
the boundary vertex changes the number of electron lines by Ag = 2. Therefore, in addition 
to the end point configurations Fig.0(a,b) at work for the infinite system, the configurations 
Fig.|2](c,d) are also allowed. 

The diagram Fig.|7], in combination with the bare Green function, determines the density 
of states in the pure system: 



Po(x,e) = po \ 1 - cos 



2[k F + -)x 

v 



(56) 



The 2fc F -oscillation of the density of states is due to the fact that the single-electron eigen- 
functions are standing waves and is not of great interest. We can get rid of this oscillation 
by adopting the following averaging procedure 

1 f 

p(x, e) -> p A (x, e) = — I dyp(x + y, e) . (57) 

-A 

Such averaging over the spatial coordinate has also been discussed by Aplil. The parameter 
A in Eq.(|57|) is chosen in such a way that 

-!- < A < min (l, -) . (58) 



(59) 



With this choice of A the 2kp oscillation in fl56|) averages away so that 



pA{x,e) = po 



O 



1 



k F A 



It is important to notice that all the diagrams with the end point configurations 
Fig||(c,d), however complicated, have an additional rapidly oscillating factor exp(±2kpx) 
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as compared to the diagrams with end point configurations Fig.[|(a,b). Upon averaging 
the former diagrams acquire [as (p9|) l a small factor of the order of ~ l/(kpA). We shall 
therefore neglect these diagrams. (We also observe that the remaining diagrams do not 
depend on A, so we drop this subscript in what follows.) 

Defining the Berezinskii blocks in full analogy to the translationally invariant case, we 
rewrite the expression (^) for the density of states in the form 

oo 

p(x, e)/p = 1 + 2Re ]T R m (x)R m (x) . (60) 

771=1 

[The mirror imaged diagrams like Fig.|5| remain such even if the boundary vertex is involved 
(with respect to the horizontal axes). Hence the factor of 2 in (pD|).] 

The right-hand-side blocks R m (x) are not influenced by the boundary at all. Hence they 
are identical to those for the translationally invariant system (Section |TV| ). As to the left- 
hand-side blocks R m (x), the presence of the boundary vertex on the left can not alter the 
process of constructing this diagrams by adding impurity vertices from the right. Hence -R m 's 
satisfy the same differential recurrence relation as they do for the infinite system: Eq. (p4|) . 
What is new in the open boundary case is that a boundary condition at x = comes into 
play. Indeed, consider the quantity R m (x —> 0). Tending x to zero means squeezing the 
diagrams R m (x) so that only 2m boundary vertices (and no impurity vertices) contribute 
to R m (0). Hence the boundary condition 

R m ( x = 0) = 1 for all m . (61) 



Next we define [c.f. p7 | 



R m (x) = (-l) m e^ mx R m , 

R m {x) = {-l) m e-^ mx L m {y) . 1 ] 



The Berezinskii blocks i? m are defined in Section [TV] [Eq.(|36|)1, while the blocks L m (y) satisfy 

= m {isL m - [AmL m - (2m - l)L m _i - (2m + l)L m+1 )} , ^ 
L m (y = 0) = 1 , 



Equation (p0|) now reads 



Po 



l + 2ReY,L m (y)Rm- (64) 



777=1 



Equation (|63| ) [in combination with Eq.(|6~4|) and Eq. (|36|) 1 determines, in principle, the 
exact density of states. Unfortunately, we failed to obtain the solution to the differential 
recurrence relation (|63|) in a closed form. [Upon making the Laplace transform with respect 
to y, one can reduce (|x|) to an algebraic recurrence relation. That is, however, not any 
simpler than the relations (0,^) for Q^ a ' b \] On the other hand, we are primarily interested 
in the low-energy (s C 1) behavior of the density of states. So, as in Section [TV], it is 
sufficient to consider the limit (0). For large m, the relation fl6"3|) reduces to the differential 
equation 
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2dL _ d ( dV 

z dy dz \ dz 



- zL , (65) 



(We recall that z = \/—2ims.) Since the operator on the right-hand-side of (|65"D is the 
operator defining the modified Bessel equation, we again (as in Section [TVp use the Lebedev- 
Kontorovich transformationil. It is convenient to define 



oo 

F(y,T) = J^K iT (z)L(y,z), 

oo (66) 
L{y,z) = h I dTTsinh(7TT)K iT (z)F(y,T) . 





Applying the transformation (p6|) to the equation (pq) we find that the function F(y, r) 
satisfies the equation 

dF(y,r) 1 2 

-^y— = -2 TF ^ T) - (6?) 
The solution of this equation is 

F(y,r) = F (r)e-y^/ 2 . (68) 
The function Fq(t) = F(0,t) is determined by the boundary condition ( |63D 

Fo(r) = f -K ir (z)L(0, z) = - , ^ . (69) 

J g 2r smh(7rT/2) 

Combining ( |69| , |68D , and (^) we find the function L(y, z) in the form 

L(3/,z) = - / drcosh(7rr/2)e- yT /2 K iT (z) . (70) 

7T JO 

This equation can be re-written in a simpler way as 

FY 7 / x 2 \ 

L(y,z) — \\ — / dx cos{z sinh x) exp . (71) 

V ^ j[ V 2 2// 

Substituting this expression into the formula flB"^]), replacing the sum over m by the integral 
over q and using the asymptotic expression (|38|) for -R m 's, we find the density of states in 
the form 

p(y, 8 )= , 2 ?*J {V)- (72) 
sm (1/s) 

The function f(y) is given by 

H2~ °f dx ( x 2 \ 

^-v^/^^ (l - xtanhx)ex H"%J • (73) 

It is plotted in Fig.H. The limiting values of the function / are 



18 



1 + 0(y) , for y < 1 

/<! " H (i) I/2 + ofe 3/2 ).fo^»i <74) 

As follows from (|72|), close to the boundary the local density of states turns out to be 
logarithmically enhanced as compared to the bulk one [c.f. fl53[)] 

p{y = V,s) = ^- (75) 
sin (1/sJ 

(see also Appendix B). Away from the boundary this enhancement is exhausted, for the 
function p(y,s) decays as (y ^> 1): 



sin (1/s) ^ 

It crosses over to the bulk value (|53|) at the distances y of the order of ln 2 (l/s) ~ L s . Thus 
the boundary influences the local density of states on the scale of the order of the correlation 
length. A similar conclusion has been reached by AP on the basis of the boundary density 
of states distribution function^ (we refer to the case when the excitations can not leave 
the sample). In our case, the difference with the standard localization problem is that 
perturbations caused by the boundary decay into the bulk on the scale L s , not A s , and that 
they follow a power law, not an exponential decay. 



C. Density of states by renormalization group 

In the Introduction we have referred to real-space renormalization group (RG) calcu- 
lations for random spin models. It is worth to investigate how this alternative approach 
works for calculating the local density of states close to the boundary. We use the real 
space version of the RG, introduced by Dasgupta and Mai and extended by Fisheri, to cope 
with disordered exchange constants in an Heisenberg chain. In fact, the random Heisen- 
berg model in the anisotropic XY limit and at zero magnetization is equivalent to a model 
of spinless fermions with random nearest-neighbor hopping integrals at half-filling. As we 
said, the random hopping has a particle-hole symmetry at half-filling, which prevents the 
generation of forward scattering processes. Therefore, this model is exactly equivalent to 
the random backscattering model. Moreover, it was shownS that the spin anisotropy is not 
a relevant parameter, and that in the whole range < J z < J x = J y the physical behavior 
is unchanged, including also the exponents of the power-law decaying correlation functions. 
Hence, we are going to consider the spin-isotropic case. 

The RG procedure consists in successive eliminations of pairs of spins, which are coupled 
more strongly than the others. The cut-off energy Q, which is rescaled downwards, is simply 
the largest exchange coupling Q = Max( J) at a given stage of the RG process. Once a pair 
of spins is decimated, an effective exchange is generated between the two spins adjacent to 
the decimated bond. The probability Q(y, T) that a bond of length y will be decimated at 
the scale fl = Q e~ r has been calculated in Ref. f| (f2 being the cut-off energy at the start 
of the RG process). 
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Let us consider a semi-infinite chain. The probability n(x, T) that the spin at site 
x G [0, oo] is still free at the scale Q satisfies the differential equation: 



dlnn(x, T) 
df 



/ dyQ(y,T)- / dyQ(y,T). 
Jo Jo 



(77) 



The two integrals on the right-hand-side imply that the spin at x gets bound to a spin at 
its left or its right. From the RG equations, one finds that Q(y, T) at T ^> 1 is given by the 
inverse Laplace transformO 



Q(y,r) 



c + io ° dz 



z 



7T 



,-Mvi lixi sinh (i/ir) 

2 oo 

* / -, \mi1 

m e 



r 3 



E( 

m=l 



1 Nm+l TT1 2 -7r 2 m 2 3 //(2r 2 ) 



(78) 



From this equation one derives that 



and, for 1 < x < T 2 , 



/°° dyQ(y,r) = i 
jo 1 

^dr' f dyQ(y,T') * ~m 
l jo 2 



The limiting expressions of n(x, T) can easily be obtained from the above equations: 



n(x, r) 



for 



x > r z 



i 



for i < x < r 2 



for 



x < 1 



In the fermion language, n(x, T) is proportional to the local density of states p(x, e) integrated 
up to the energy e = fi e ~ r - We finally obtain 



p(x,e) 



|e|ln d M/|e 



^|6|ln 2 (l/|e|) 



for x ^> L e 
for 1 < x < L f 



I lelhWl, 



for 



x 



where the cross-over length L e = 21n 2 (e)/7r 2 . The resulting density of states is therefore 
in perfect agreement with the previous calculation. An important observation is that from 
the RG analysis it comes out quite naturally that the relevant correlation length in this 
problem is L t . This is one of the main achievements of this approach, as we are going to 
discuss in the following section. Before concluding this section, it is worthwhile to briefly 
discuss the enhancement of the local density of states close to a boundary in the light of the 
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RG approach. The RG provides a quite simple picture of the ground state of the random 
Heisenberg model, in terms of the so-called 'random singlet state's. With this we intend a 
state where each spin is coupled into a singlet with another spin, and the longer the distance 
between the two coupled spins the lower the excitation energy of the singlet. Hence, it is 
more or less obvious that a spin close to a boundary has less probability to get coupled to 
another spin, which explains the enhancement of the density of states at low energy. 



VI. DISCUSSION AND CONCLUSIONS 

To summarize, in this paper we have investigated the single-particle Green function and 
the local density of states in the presence of boundaries for the random mass Dirac problem. 
Let us discuss these findings in reverse. 

For the density of states of the random mass Dirac fermions is strongly energy dependent 
(Dyson singularity), it is not a surprise that the presence of the boundary is essential. It is 
perhaps less obvious that the density of states close to the boundary is strongly enhanced 
as compared to the bulk onei! 

p{x = 0,s)/p(s)~ iln(l/s). (79) 

Regarding the spin systems, this result is relevant for the situation when one considers 
the effect of especially strong impurities (like lattice defects) in an already doped sample. 
Indeed, the Dyson singularity itself readily yields peculiar thermodynamic properties^'". 
For example, the linear magnetic susceptibility diverges as xCO ~ 1/Tln 2 (l/T)t3. Heeding 
( [79]) , one expects (larger) contributions of the order of 

X (T) ~ l/7Tn(l/T) 

to the magnetic susceptibility due to strong impurities and, of course, sample boundaries 
(the latter are especially important for powder samples on which most of the experiments 
are presenty being conducted!!!) . 

According to the result (0), the single-particle Green function behaves as 

n / 1 \ 3 / 2 

G + (x,s)~ 2 P ° , - -) , (80) 



s In (—is) \xj 

in the region i <C i <C i s . (In the region x ^> L s we expect the Green function to cross 
over to a standard exponential decay.) It is worth noting that the Green function (j30l) has 
a large prefactor. Therefore it is actually larger than the bare Green function G*( 0+ ) ~ p 
on the whole scale 1 C i < L s . The conclusion is that not only the density of states is 
large (that is not surprising for the zero-modes are packed around zero energy), but also the 
correlations are strong. In fact, the one-particle correlations in the disordered system are 
stronger that those in a pure system (without the gap). This is only true below a certain 
energy scale, which vanishes with the impurity concentration. 

The length scale L s ~ I In 2 s, below which Q8"H| ) is valid, can be inferred from the analysis 
of the boundary effects. The same scaling relation of the length x versus the energy s, 
x ~ I ln 2 (s), was found in Ref. [12] by means of a replica method analysis of the Dirac fermion 
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model, and also was derived by the RG investigation of the random Heisenberg model in 
Ref. £|. Surprisingly, it differs from the scaling relation that the energy dependence of the 
localization length would suggest, namely x ~ A(s) = l\ lns|. This is an evident signature 
that the typical behavior differs from the average behavior, the difference becoming large as 
s — > 0. In other words, some physical quantities are dominated by very rare events. 

Notice that L s does not come out in any simple way from our asymptotic calculation 
of the Green function. The likely reason is that, once we solve the recurrence relations for 
the Berezinskii blocks under the restriction fl37|), we are implicitly assuming that L s ^> x. 
This denies the access to the region x ^> L s where the Green function decays exponentially. 
This issue ought to be further studiecH. In fact, in the case when we knew the x L s 
asymptotics, that is for the density of states in the presence of a boundary, we have been 
able to identify L s . 

As we mentioned in the Introduction, also the RG method gives access to the correla- 
tion functions, although at equal times. In particular, the equal-time spin-spin correlation 
function is known to decay at large distances as! 

(S\x)S\0)) ~ 1 
x A 

{% = x,y,z). One can showH that also the equal-time Green function in the equivalent 
random hopping tight-binding model decays as 1/x 2 . This result is, in fact, compatible with 
(jSOf). Indeed, (|80|) is only valid for i C I s , or, equivalently, for s C s(i) = exp(— aJ~x/1). 
Therefore, neglecting exponentially decaying terms, we find 

r s ( x ) r>n r s ( x ) 1 

G + (* = °,*)~i ««?♦<«,.) ~^jf -^ r - 2 , (81) 

in agreement with the RG results. Hence our results not only reproduce the equal-time 
behavior of the correlation functions found by RG, but also allow us to determine their 
energy dependence. 

According to GM0, the two-particle correlations follow the same pattern as (|80|). The 
fact that the disorder actually enhances correlations (at low energy scale) throws the light 
onto the experimental findings of an antiferromagneting ordering in the spin-Peierls com- 
pounds building up upon dopingS. Indeed, the system without the spin-phonon interactions 
is known to be unstable with respect to the antiferromagnetic ordering. In a pure system, 
however, the latter is prevented by the spin gap. The doping effectively creates states in 
the gap and induces correlations that are even stronger than those in a pure system without 
the spin gap. Hence the antiferromagnetism is promoted by doping. This matter, though, 
is to be investigated in more detail, especially with respect to the energy scales involved 
(we expect the relevant scale be associated with the fluctuations of the density-density 
correlator )0'0. 

To conclude, in this paper we have analyzed single-particle properties of the random 
mass Dirac problem in some detail. The lack of knowledge of the distribution functions (for 
the density of states, etc.), however, calls upon further exploration of the problem. 
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VIII. APPENDIX A 

In order to evaluate the integrals in ([18|, |49|) , which read 

oo oo 

F{y) = J dzK (z) J drrtanh (^j K tT (z)e~^ y , (82) 
o o 

we make use of Nicholson's integral representation of the product of MacDonald functions!! 

oo 

K (z)K iT (z) = 2 J dtcos(Tt)K iT (2zcosht) . (83) 



Upon substituting fl83|) into fl82| ) one observes that the ^-integral and the t-integral, if taken 
in succession, are the table onesil. Hence (|5l|) . 

IX. APPENDIX B 

Exactly at the sample boundary the local density of states can be explicitly evaluated. 
Indeed, owing to the boundary condition (pTj) , the formula ( ]oTi| ) for the density of states 
simplifies so that the summation over the number of lines is easily performed [c.f. (|36|)1: 



p(x = 0, e)/p = 1 + 2Re £ R m = -1 + 2Re <j Kq{ _ w/a) J dt \j -^e^ 2 



m=l 



The above integral is a table oneE3, so that we finally obtain 



p(x = 0,e)/ Po = Re{ ^ ^ 1 (84) 



K (-is/4) 

[the s — > limit of this formula is in agreement with Eq.(|72|)l 
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FIGURES 

FIG. 1. The elementary vertices of the Berezinskii diagram technique contributing the factors 
(a) -1/(2/), (b) -l/l, (c) -(l/l)e W [-i(e/v)x], (d) -(1/Z) exp [+i(e/v)x]. 

FIG. 2. The configurations of the electron line end points for the single-particle Green function. 

FIG. 3. This figure illustrates the definition of Berezinskii blocks for the diagrams with end 
point configuration shown in Fig. ||(a). 

FIG. 4. The process of joining the elementary vertices for deriving the recurrence relation for 
the right-hand-side blocks R m (x). 

FIG. 5. Two fourth order mirror imaged (coinciding end points) diagrams contributing to the 
Green function. 

FIG. 6. The plot of the function F(y) determining the spatial decay of the single-particle Green 
function. 

FIG. 7. The new element of the diagram technique in the presence of the boundary - the 
boundary vertex. 

FIG. 8. An example of a diagram involving the boundary scattering. 

FIG. 9. The plot of the function f(y) determining the spatial decay of the local density of states 
away from the boundary. 
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